False Vacuum Decay Induced by Particle Collisions 
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The semiclassical formalism for numerical calculation of the rate of tunneling transitions 
induced by N particles with total energy E of order or higher than the height of the barrier 
is developed. The formalism is applied to the induced false vacuum decay in the massive 
four-dimensional — XcjA model. The decay rate, as a function of E and N, is calculated 
numerically in the range 0.4 < E/E sp h < 3.5 and 0.25 & N/N sp h Ss 1.0, where E sp h and 
N sp h are the energy and the number of particles in the analog of the sphaleron configuration. 
The results imply that the two-particle cross section of the false vacuum decay is exponentially 
suppressed at least up to energies of order X0E sp h- At E ~ E sp h, this exponential suppression 
is estimated as about 80% of the zero energy suppression. 



PACS numbers: ll.10.Jj, ll.10.Kc, 12.38.Lg 



I. INTRODUCTION 

Non-perturbative effects related to tunneling play an 
important role in many field theory models. Two widely 
known examples are the decay of the metastable (false) 
vacuum in scalar models and winding number transitions 
in sigma models and gauge theories at weak coupling. At 
low energies, these processes are well described by the 
semiclassical approximation which relies upon the exis- 
tence of classical (imaginary time) solutions interpolating 
between initial and final states. In the two examples men- 
tioned above the interpolating solutions are bounce ||] 
and instanton || , respectively. The Euclidean action of 
these solutions determines the exponential part of the 
transition rate. It is inversely proportional to the small 
coupling constant; the transition rate is thus strongly 
suppressed. In realistic situations, the initial state may 
contain particles. If these particles have low energies, 
they can be taken into account perturbativcly and merely 
change the pre-exponential factor. 

The situation is different if particles with parametri- 
cally high energy are present in the initial state. An ex- 
ample relevant to what follows is high energy scattering 
in the false vacuum or in the presence of the instanton. 
If the energy of the collision is of order 1/A, where A is 
the small coupling constant, then the application of the 
semiclassical approximation is not straightforward, de- 
spite the fact that the two-particle initial state is pertur- 
bative. The probability of induced tunneling transition 
may, and in fact does depend exponentially on the energy 
of collision. Quantitative description of this effect is the 
main purpose of this paper. 

The exponential dependence on collision energy is pre- 
cisely what one expects on physical grounds. The typical 
situation is illustrated in Figjl] which shows a generic 
profile of the potential energy as a function of field for a 
model with false vacuum decay. The false vacuum = 
is separated from the true one by the potential barrier 



whose heightS E sp h is parametrically m/X, where m is 
the mass scale in the model. The intuition based on 
quantum mechanics suggests that at energies of order 
E sp h the tunneling suppression factor may be substan- 
tially reduced or even completely absent. This quantum 
mechanical intuition, however, should be translated to 
field theory with care. 
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FIG. 1. Generic picture of false vacuum decay induced by 
particle collisions. The wavy line represents an excited state 
above the false vacuum. 

The situation is much the same in the case of instanton 
transitions, at least in models which possess a mass scale 
at the classical level. An example is the Electroweak The- 
ory where different topological sectors are separated by 



1 Throughout this paper we somewhat loosely use the term 
"sphaleron" not only for the particular static solution in the 
Electroweak Theory [3|, but for any static solution sitting 
in unstable equilibrium on top of the potential barrier and 
thus characterizing its height. In the problem of false vacuum 
decay the sphaleron is actually the critical bubble Q ; E sp h in 
the text stands for the sphaleron energy. 
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the potential barrier of the height E s ph ~ lOTeV [[|. 
It was found || that tunneling between different sec- 
tors is enhanced in the presence of high energy parti- 
cles. Whether the exponential suppression disappears 
completely at sufficiently high energy is still an open 
question. Were it possible to observe such transitions in 
collider experiments, the events would look quite spec- 
tacular due to accompanying baryon and lepton number 
violation || (for recent review on baryon number viola- 
tion see ref. 0). 

As was first noted in refs. J5j, at relatively low energies 
the corrections to the tunneling rate can be calculated 
by perturbative expansion in the background of the in- 
stanton (bounce). Further studies showed that the actual 
expansion parameter is E/E sp h [^P Jlo| , [ll|| and the total 
cross section of the induced tunneling has a remarkable 
form 



a tot (E) ~ exp ^-jF HG (E/E sph )^ 



(1.1) 



where the function! Fug ( e ) is a series in powers of e (for 
a review sec |l2"|). 

While the perturbation theory in e = E/E sp h is lim- 
ited to small e, the general form of the total cross section 
implies that there might exist a semiclassical-type proce- 
dure which would allow, at least in principle, to calculate 
Fhg(^) at e <; 1. Since the initial state of two highly en- 
ergetic particles is not semiclassical, the standard semi- 
classical procedure does not apply and a proper gener- 
alization is needed. Such generalization was proposed in 
refs. @||,|||. The formalism of refs. (l|jl|Jl| reduces 
the calculation of the exponential suppression factor to 
a certain classical boundary value problem, whose ana- 
lytical solution is not usually possibleS. The purpose of 
this paper is to develop numerical techniques for calculat- 
ing the function Fjjg(E / E sp h) and to demonstrate how 
this techniques works in a simple model. As an example 
we consider the induced false vacuum decay in the — A</> 4 
theory with the mass term (in a different context, the 
problem of induced false vacuum decay was previously 
addressed in refs. |i"^Jt6| ]). The rate of induced tunneling 
transitions should be calculable by the same method in 
more realistic models as well. 

The motivation for this calculation is the following. Al- 
though general arguments imply [Q that the cross sec- 
tion of induced tunneling is exponentially suppressed at 
all energies (our preliminary numerical results Jl9| l agree 
with this statement), in realistic models with finite cou- 
pling constant the possibility to observe induced tunnel- 
ing transitions depends on the numerical value of the sup- 
pression factor. The only way to estimate this value is by 



2 The subscript HG stands for "holy grail" p|. 

3 In some simple cases a non-trivial information about the 
function Fmg(E / E sp h) can still be obtained by analytical 
methods 



performing direct calculations. This may be of particular 
importance in condensed matter models where coupling 
constants are often not very small and even low energy 
transitions are observable. 

The semiclassical approach proposed in refs. El , , ^5| 
is based on the conjecture that, with exponential accu- 
racy, the two-particle initial state can be substituted by 
a multiparticle one provided that the number of particles 
is not parametrically large (although not proven rigor- 
ously, this conjecture was checked explicitly in several 
orders of perturbation theory in e [ p0[ ) . The few-particle 
initial state, in turn, can be considered as a limiting case 
of truly multiparticle one with the number of particles 
N = v/\ when the parameter v is sent to zero. For the 
multiparticle initial state the transition rate is explicitly 
semiclassical and has the form 



a(E, N) ~ exp<^ ——F(e, v) 



A 



According to the above conjecture, the function Fhg{z) 
can be reproduced in the limit v — > 0, 



lim F(e, v) 



F HG (e). 



Therefore, although indirectly, the function Fhg( 6 ) is 
also calculable semiclassically. 

Within the conventional semiclassical framework, the 
function F(e, v) is determined by the action evaluated on 
a particular solution to the classical field equations with 
certain boundary conditions Jl5| . In this formulation, 
the problem allows for numerical study. Namely, one 
can solve the corresponding boundary value problem nu- 
merically and calculate the function F(e,v), which then 
can be used to extract information about Fhg{^)- In the 
present paper we follow this general strategy. On the lat- 
tice, it is not possible to actually reach the point v = 0, 
as the solution to the boundary value problem becomes 
singular in this limit. It is important, however, that the 
function F{e 1 v) is monotonically decreasing function of 
v at fixed e |l3|,[l4]], so that the point v = corresponds 
to maximum suppression at given e. If one finds that 
F(e, v) > for some z/, then 



F 



H G (e) > F(e,u) > 



and the two-particle cross section is exponentially sup- 
pressed. The energy E* as which the two-particle cross 
section would become unsuppressed is the energy at 
which the line F(e, v) = would cross the v = axis. 

In the — A</> 4 model with the mass term (described in 
detail in Sect. 3), we calculate the function F{e 1 v) numer- 
ically in the range 0.25 < v < 1 and 0.4 < e < 3.5, where 
v and e are normalized to their values at the sphaleron. 
We perform extrapolation to v = and estimate the 
value of the function Fhg(() at e ~ 1, i.e. at energies 
of the order of the sphaleron energy. We find that at 
these energies the zero energy exponential suppression of 
the two-particle cross section is reduced by only about 
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20%. By extrapolation of F(e, v) to higher energies we 
set lower bound for the energy E* at which the func- 
tion Fhg(^) may become zero and thus the two-particle 
cross section may become unsuppressed. We find that 
E* > 10 E sp h- Hence, the tunneling transitions induced 
by particle collisions are exponentially suppressed well 
above the sphaleron energy, at least in the model we con- 
sider. 

The paper is organized as follows. In Sect. II we discuss 
in more detail the definition of the multiparticle proba- 
bility cr(E, N) and the use of the semiclassical approxi- 
mation for the calculation of the corresponding exponent 
F(e,u). We rewrite the formalism of ref. in a form 
suitable for numerical calculations. In Sect. Ill the — A</> 4 
model with the mass term, as well as general features 
of the false vacuum decay in this model, are described. 
We then turn in Sect. IV to analytical calculation of the 
function F(e, v) at e <C 1 and at those values of v which 
maximize the transition rate. This is done for comparison 
with numerical results. In Sect.V we present the results 
of numerical calculation of the function F(e, v). Sect. VI 
contains discussion and concluding remarks. The details 
of numerical techniques are given in Appendices A and 
B. 



II. SEMICLASSICAL CALCULATION OF 
INCLUSIVE MULTIPARTICLE PROBABILITY 

In this Section we review the the formalism of refs. 
^3| , ^4|Jl5| and rewrite it in the form suitable for numerical 
calculations. 

The inclusive multiparticle probability a(E,N) dis- 
cussed in the Introduction is defined as follows, 



a(E,N)=Y / \(f\SP E pN\t) 



(2.1) 



where S is the S-matrix, Pe,n are projectors onto sub- 
spaces of fixed energy E and fixed number of particles TV, 
respectively, while the states \i) and |/) are perturbative 
excitations above two vacua lying on different sides of the 
barrier. a(E, N) can be interpreted as the total proba- 
bility of tunneling from a state of energy E and number 
of particles N, summed over all such states with equal 
weight. 

The advantage of dealing with a(E, N) instead of 
atot(E) is that the r.h.s. of eq.( |2.l[ ) can be written in 
the functional integral form, in which the semiclassi- 
cal approximation is equivalent to the saddle-point in- 
tegration. The double path integral representation for 
a(E,N) reads p| 



<r(E,N) = j d6dTda k dakdbkdb^d<j>(x)d<f>'(x)exp{-iN6 

—iET - J dka k ale- %9 - ZbJkT - J dkb k bt 
fa) + B f (bl cj> f ) + B? (oik. 



+B}(b. k , <f>' f ) + iSty) - (2.2) 
where the boundary terms Bi and Bf are 



Bi(a k ,<pi) = i Jdk 



-Wk<M*0<M-fc) - aka-ke" 



-2iuj k Ti 



Bf(bi,<t> f ) 



+2^/ 7 M ^ e- lulkT >a^ ^ (j) i {k) 
1 



dk 



^4>m<t> f {-k)-bib*_^^ 



+2 % /2^ k V WfeT '6 k / (-fc) 



(2.3) 



Here <f>ij(k) are the spatial Fourier transforms of the 
field at initial and final times Tj and Tf, respectively. 
The limit Tij — > =poo is assumed at the end of the cal- 
culation. The complex integration variables a k and 6£ 
come from coherent state representation of initial and 
final states |2lj] ; they are classical counterparts of annihi- 
lation and creation operators. The integration over these 
variables substit utes the summation over initial and fi- 
nal states in eq. fl2.ip . The functional integrals over <f>(x) 
and 4>'(x) come from the amplitude and complex con- 
jugate amplitude, respectively. The integrations include 
the boundary values tfiij and <f>\ f. 

From eq.( [2.2| ) it follows immediately that in the wea k 
coupling limit a(E,N) has the semiclassical form fl2.l| ). 
Indeed, changing the integration variables 

{a, b, (j), (j)'} -> ^{a, b , 0'} 

and taking into account that after this transformation the 
action becomes proportional to 1/A, we arrive at eq.( [2.1| ) 
where e = XE, v = \N and the function F(e, v) is de- 
termined by the saddle-point value of the exponent in 

Let us now discuss the saddle-point equations for the 
integral (2.2). We will see that these equations reduce 
to certain boundary value problem for the fields <p and 
if)'. The variables a k , a£, 6 k and do not play any 
role in what follows. Moreover, they enter the exponent 
quadratically. Integrating them out we find 



a(E,N)= / d6dTd(j>{x)d<i>' (x) JJ<5(^/(fc) -</>/(&)) 



x exp<^ -iN6 - iET + iS(<j>) ~ iS((j>') 



where 



\ J *73^(( 1 + 7k)[<Mfc)M-fc) 

~ 4 7k ^(fe)^(-fe)) V (2.4) 



7k = e ie+i ^ T . 



3 



The important feature of the representation (2.4) is that 
the exponent in the r.h.s. contains only the action and 
the boundary values of the fields. Thus, the discretization 
of this expression is straightforward. 

Let us turn to the saddle point equations. Varying the 
exponent with respect to the fields <f>(x) and <p'(x) we find 



SS SS 



0, 



(2.5) 



i.e. the usual field equations. The boundary conditions 
for these equations come from the variation with respect 
to the boundary values of the fields. At t = Tf, because 
of the ^-function, the variations are subject to the con- 
straint 8(f) f {x) = 6<p' f (x). Since SS/S(j)(T f ,x) = (j){T } ,x) 
we obtain 



j>(T f ,x)=j>'(T f ,x) 
cf>(T f ,x)=4>'(T f ,x) 



(2.6) 



Thus, in the final asymptotic region the saddle-point 
fields 4> an d <j>' coincide. Note that the dependence on 
Tf cancels out in the difference S(4>) — S(4>'). 

The variation with respect to 4>i and </>^ leads to two 
equations which can be written in the following form, 



i(j)i{k) +w k 0j(fc) = 7k + ^k^( fe ) y 

-i0i(fc)+«k&(fc) = — (-ii>'i(k) +w k ^(fe)) • (2.7) 

7k ^ / 

Let us check that these boundary conditions imply the 
independence of the exponent in cq.(2.4) of the initial 
normalization time Tj. The action depends on T t ex- 
plicitly, while the boundary term in cq.( |2.4| ) depends on 
Ti through the boundary values of the fields, <f>i and <\)\. 
Assuming that the fields linearize and satisfy free field 
equation in the initial asymptotic region, and integrating 
by parts in the actio n, t he li-dependent contributions in 
the exponent in cq.(2.4) read 



dk 



[(i + T£)(&(fc)&(-fe) + #(fcM(-fe)) 



-4 7k ^(fe)^(-fe)] + ifc(fc)&(-*) - i#(fc)#(-fc) 



By making use of eqs.(2.7) it is straightforward to see 
that the integrand in this expression vanishes. 

Finally, there are two saddle-point equatio ns which 
come from the variation of the exponent in eq^^) with 
respect to 9 and T. They read 



E 
N 



i (fc)-7k^(fc)][^(-fc)-7k0 l (-fc)], 



2o;^7krffc . 

(i-7k0 2 ' 

2c,k7kdfc [4(fc) _ ^ mm _ k) _ 7k ^ i( _fe)]. 

(2.8) 



These equations determine the saddle-point values of 9 
and T in terms of E and N. 



The initial boundary conditions (2.7) simplify when 
written in terms of frequency components. In the initial 
asymptotic region where <p and </>' are free fields, we can 
write 



4>{x) 



dk 



= {/ k e 



--icjfci+zkx i * iwkt— ikx 



<P ( X ) = I 777. \tn \/k P 



5ke 

/* „. 
3ke 



}■ 



(2.9) 



Eqs.( [2.7| ) then become 



/k = 7k / k 
* 1 

7k 



k- 



(2.10) 



while the saddle-point equations (2^8) read 
E = I dku) k f k g^, 



N = J dkhgl 



(2.11) 



One may recognize the usual expressions for the energy 
and the number of particles contained in the free classical 
field, /kSk being the occupation number in the mode with 
spatial momentum k. 

Since, as follows from eq.(2.6), the fields <j> and (/)' co- 



incide in the final asymptotic region and thus coincide 
everywhe re, the re a re solutions to the boundary value 
problem (2.5) - (2/7) which have 7k = 1 and 4> = 4>' in the 
initial asymptotic region. These solutions do not describe 
tunneling transitions — for them the exponent in eq. ( |2.4| ) 
vanishes. They correspond to classical propagation over 
the barrier (recall that we require that the initial and fi- 
nal states are on different sides of the barrier) and clearly 
exist only at E > E sp h ■ If such over-barrier solutions ex- 
ist for some values of E and N, then the corresponding 
multiparticlc probability is definitely not exponentially 
suppressed. Conversely, if the two-particle cross section 
are not exponentially suppressed at some energy, then at 
this energy there should exist over-barrier solutions with 
arbitrarily small N. This observation provides a strategy 
of search for unsuppressed induced transitions p2j . Our 
formalism is not suitable for stu dying th ese solutions as 
the boundary conditions eq.(2.6) and eq.(2.7) are degen- 
erate in the case (f>i = (fy, 7k = 1. 

The solution to the boundary value problem (pT^)-(p?7|) 
which describes tunneling transition is obtained when, in 
the initial asymptotic region, (j) and <j)' are treated as val- 
ues of the same field <f> taken on different sheets in the 
complex time plane. The existence of such solutions is 
suggested by the analogy with barrier penetration in one- 
dimensional quantum mechanics (below we show by nu- 
merical calculations that analogous solutions exist in the 
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field theory as well, and argue that their interpretation 
as tunneling solutions is consistent). 



per period. Thus, we recover the standard WKB transi- 
tion probability exp{i;T - S E (T)}. 




FIG. 2. Barrier penetration in one-dimensional quantum 
mechanics. The classical solution with given energy, consid- 
ered in the complex time plane, describes both evolution in 
the classically forbidden (BC part) and classically allowed re- 
gions (AB and CD parts). 
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FIG. 3. The contour in the complex time plane where the 
boundary value problem (6)-(8) is formulated. Crossed circles 
represent singularities of the field. 



To illustrate this idea consider barrier penetration in 
one-dimensional quantum mechanics with the potential 



V(q) = 



1 



ch 2 <7 



shown in Fig.g. At energy E < 1 there exist two clas- 
sically allowed regions, AB and CD, separated by the 
classically forbidden one, BC. The classical solution with 
energy E is given by the following equation, 



E 



1 — E 



shq = ch(V2Et). 



(2.12) 



As defined by this equation, q(t) is an analytic (except 
for isolated points, see below) function which is periodic 
in Euclidean time with the period T — 2n j \J2E and real 
at purely imaginary time Rei = where it describes os- 
cillations between the two turning points B and C. It 
is also real on the lines Imi = nT/2, n 6 Z where it 
represents the motion in the classically allowed regions: 
in the region CD for even n and in the region AB for 
odd n. The entire evolution from q = — oo to q = +oo, 
including tunneling, corresponds to the contour ABCD 
in the comp lex time plane as shown in Fig|| The so- 
lution ( 2.12 ) taken on the contour ABCD contributes 
to the amplitude, while taken on the complex conjugate 
contour A'B'CD it contributes to the complex conjugate 
amplitude. These two parts of the solution are quantum- 
mechanical analogs of the fields 4> and <p' '. In perfect anal- 
ogy to eq. (|]4|) , the transition probability can be written 
as exp{iSABCD — ^Sa'B'cd + ET}, where subscripts in- 
dicate the contour along which the action is evaluated. 
In the sum of the actions the Minkowskian parts of the 
contours cancel, while the Euclidean parts add and give 
iSsB' = —Se(T), where Se(T) is the Euclidean action 



For this picture to be consistent, there must be sin- 
gularities in the complex time plane which prevent the 
deformation of the contour ABCD to the real time axis. 
In particular, the initial asymptotic region A should be 
separated from the real time axis by a cu t. It is straight- 
forward to check that the solution (2.12) indeed has the 
right structure of singularities. The function q(t) itself is 
finite, but the derivative 



dq 



is infinite at the points t* satisfying 
E 



1 - E 
i.e. at the points 



ch z V2EU = 0, 



(2.13) 



(2.14) 



Note that, as follows from eqs.( 2.12| ), ( [2.13] ), the corre- 
sponding values satisfy 

dig* = 0. 

Th e typ e of the singularity can be found by solving 
eq.( 2.12j ) in the vicinity of t he sin gular point, t = t* +St, 
q = q« + Sq. Expanding eq.(2.12) to the leading order we 
find 

(Sq) 2 = ±2V26t. 
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Thus, the singularities are the square root branching 
points. The cuts relevant for our boundary value problem 
are shown in Figj|[ 

In the spirit of this quantum mechanical example, we 
interpret the fields 4> and <j)' in the initial asymptotic 
region as values of the same field <fi taken on different 
sheets in the complex time plane. The fields 4> and <£>' are 
obtained by analytical continuation (in the initial asymp- 
totic region) to the real time axis from the two complex 
conjugate contours ABCD and A'B'CD. It is convenient 
to reformulate the boundary value problem directly in 
terms of the fields on these contours (note that the ana- 
lytical continuation in the initial a sym ptotic region can 
be done explicitly by means of eqs.(2.9)). Let us assume 
for a moment that the saddle point values of T and 8 are 
purely imaginary and positive. Then 



7k 



(2.15) 



and the analytical continuation in eqs.(2.7)-(2.11) from 



the real time axis, where they are originally formulated, 
to the contours ABCD and A'B'CD results in the substi- 
tution of 7k by 



7 



In other words, eqs.(2.7)- (|2.11 ) remain valid if simulta- 
neously with the analytical continuation the substitution 
7k — > 7 is performed. Note that in this formulation the 
parameter T enters the boundary conditions implicitly 
as the total amount of Euclidean evolution. It is worth 
noting also that, except for the asymptotic region, the 
contours ABC and A'B'C can be arbitrarily deformed, 
provided the singularities are not crossed. 

In gen eral , there may be several saddle points in the 
integral (2.4). In this case, the physically relevant one 
is that continuously connected to the saddle point which 
emerges in the low energy perturbation theory |l4|] . The 
perturbative saddle point has purely imaginary 8 and T. 
Moreover, it obeys the property 



(2.16) 



In particular, <p is rea l on the- real time axis. Unless a 
bifurcation occurs, these properties should hold for the 
non-perturbative saddle point as well. In numericaLcal- 
culations we have not found any bifurcation pointsu, so 
in what follows we assume that these are absent and con- 
sider saddle points obeying eq.( [2.16 ). 

Eq.(2.16) implies that on the contour A'B'CD the field 
is complex conjugate of that on the contour ABCD. Thus, 
the field is real on the part CD of the contour, and only 
the part ABC of the contour needs to be considered. In 



4 The numerical algorithm we use (see Sect. 5) stops to con- 
verge in the vicinity of bifurcation points. 



terms of the field on the contour ABC, the multiparticle 
probability reads 

a(E,N) ~expl.~F(e,v) 



—F(e, v) = NO + ET + 2 Re[iSU bc (</>)] 
A 



1-7 
1 + 7 
1 + 7 
1^7 



dfca^Re0i(fc)Re<k(-k) (2.17) 
dkujk Im<^j(fc) Im0j(— k), 



where 7 is given by eq.(2.15). The field 4> satisfies the 
following boundary value problem formulated on ABC, 



^ = 
5cf> ' 

Im0(O,x) = lm^(0,x) = 0, 

/k = e g k , 



(2.18a) 

(2.18b) 
(2.18c) 



where /k and g^ are frequency components of the field in 
the asymptotic region A along AB. Note that eq.( 2.18b| ) 
implies reality of the field on the real time axis. It consists 
of two real equations. Eq.(2.18c) also consists of two 
real conditions imposed on the field and its derivative 
in the asymptotic region A. In total, there are two real 
(one complex) boundary conditions at each end of the 
contour, so the bounda ry v alue problem is completely 
specified. Finally, eqs.( [2.1l| ) which determine the two 
auxiliary parameters T and 9, become 



N 



2 7 



, ^1— Re0 l (fe)Re0. i (-fe) 
V(l + l) 



(1-7) 2 



E = 



dkuj 2 



Im<^(fe)lm0j(-fc) ), 
2 7 



M + 7 )2 Re<fc(fc)Re<fc(-fc) 



27 



(1"7) : 



Im 6i(k) Imc 



-fe) 



(2.19) 



This is the boundary value problem we solve numerically 
in the present paper. 

The int erpret ation of solutions to the boundary value 
problem (2.18a) is the following. On the part CD of the 
contour, the saddle-point field is real; it describes the 
evolution of the system after tunnelin g. On the contrary, 
it follows from boundary conditions ( 2.18c ) that in the 
initial asymptotic region the saddle-point field is complex 
provided that 8 =/= 0. Thus, the initial state which maxi- 
mizes the probability ( |2.2| ) is not described by a classical 
field, i.e. this stage of the evolution is essentially quantum 
even at N ~ 1/A. 

The case 8 = is exceptional. In this case, the bound- 
ary condition (2.18c) reduce to the reality conditions 
imposed at Imt = T/2. The solution to the resulting 
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boundary value problem is given by periodic instanton 
of ref . ]23[ . Periodic instanton is a real periodic solution 
to the Euclidean field equations with period T and two 
turning points at t — and t — iT/2. Being analytically 
continued to the Minkowskian domain through the turn- 
ing points, periodic instanton is real at the lines Imi = 
and Imt = T/ 2 and therefore satisfies the boundary 
value problem ( [2. 18a ) with 8 — 0. Periodic instanton 
is the close analog of the quantum-mechanical solution 
discussed above. 



III. THE MODEL 

We perform numerical calculations in the model of one 
real scalar field 0, defined by the Minkowskian action 



S = 



d i x(±d^<f>-^m 2 <j ) 2 + ±\<f> i ), 



(3.1) 



where A is a positive constant. In the case m = and 
infinite volume, the state cj> = is stable classically but 
meta-stable with respect to quantum fluctuations. The 
decay of this state is described by the well-known instan- 
ton solution I2IJ2I, 



P 



A x 2 + p 2 



Due to conformal invariance, the instanton can have ar- 
bitrary size p. The instanton action, 



Si, 



8^2 
3A ' 



(3.3) 



is independent of p. At the same time, the energy bar- 
rier between the state </> = and the instability region 
depends on the configuration size and tends to zero as 
p — ► 00 . 

At m ^ 0, the conformal invariance is softly broken 
and regular Euclidean solutions with finite action are for- 
bidden by scaling argument. The action of the instanton- 
like configuration depends on its size and is minimal at 
p = 0. The decay of the state (j> = is dominated by a 
set of approximate solutions, constrained instantons p6[ , 
minimizing the Euclidean action under the constraint 
that their size is p. At p 2 m 2 <C 1, the field of the con- 
strained instanton coincides with the massless solution 



( p.2[ ) at x <C to and falls off exponentially at x <; m . 
At low energies, contribution of small size constrained 
instantons is dominant and the probability of the decay 
is determined by the factor exp(— Si ns t). 

At m ^ 0, the barrier separating the state 4> = from 
instability region (f> m/VA is finite. There exists the 
sphaleron solution characterizing the barrier height. It is 
a static 0(3)-symmetric configuration 4'sphij') satisfying 
the equation 



- -4>' + mcf> - X(j) 3 = 0, 
r 



which can be solved numerically by shooting method (see, 
e.g., ref. |27]]). The energy of the sphaleron equals 



E sp h = 18.9—. 



Since the sphaleron corresponds to the top of the poten- 
tial barrier, it is unstable and has one negative mode. 
In the model ( |3.1| ) the negative eigenvalue can be found 
numerically, 

oj 2 _ = -15.3to 2 . 

Being slightly perturbed in the direction of the negative 
mode, the sphaleron rolls down to the metastable vac- 
uum and decays into particles (i.e., becomes a collection 
of plane waves). In this way the number of particles 
contained in the sphaleron is defined p8[ . In numer- 
ical sim ul ation s it can be me asur ed by making use of 
eqs. (2.9), ( 2.11 ). In the model (3J) we found 



N 



sph 



10.5 



The peculiar feature of the model (3.1) is the absence 
of a stable vacuum state. However, this is not essential 
in semiclassical calculations. One may imagine that the 
(3.2) term ct(f> 6 is added to the potential with a small coeffi cient 
a. Such term would create a true vacuum at 4> ~ \J A/a 
but would not change equations of motion in the region 
<f> < to/ VA. In numerical calculations it is more conve- 
nient not to add such a term since the singularity in the 
final state (i.e. on the part CD of the contour of Fig.||) is 
a clear sign of the false vacuum decay. 

Before turning to the calculation of the function F(e, v) 
in this model, it is useful to note that neither of the 



two parameters entering the action (B.l) appears in the 
classical equations of motion. Indeed, after changing the 
variables according to 



m 1 x, 



7T' 



the action (3^) takes the following form, 

, 1 ,0 1 



5' 



d 4 x 



It is also clear from this equation that the coupling con- 
stant A controls the semiclassical approximation. 



IV. PERTURBATION THEORY ATe«1 AND 
9 = 

At low energies the function F(e, v) can be calculated 
analytically. For the purpose of comparison with numer- 
ical results, we perform this calculation in the case 6 = 0, 
i.e. for the periodic instanton. 
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FIG. 4. Low energy periodic instanton is a sum of con- 
strained instantons sitting at the distance T along the imag- 
inary time axis. 

At e <C 1, the periodic instanton can be approximated 
by an infinite chain of zero energy constrained instantons 
of the size p sitting at the distance T along the imagi- 
nary time axis, as shown in Fig.^. The action per period 
for this field configuration, S p , is a function of p and T. 
Since, by reflection symmetry, the point t — iT/2 is a 
turning point, <p(iT/2,x) = 0, the field is real on the 
line AB of the contour of Fi g.|I|. Thus, it satisfies the 
boundary condition eq. fl2.18c ) with 6 = 0. One also has 



2Re[iS ABC (<f>)] = S p (p,T) ) 



while the boundary term in the exponent in eq.(2.17) 
vanishes. Therefore, for the case of periodic instanton, 
the exponent in eq.(2.17) reads 



~F(e,v p (e)) = ET - S p (p,T), 



(4.1) 



where the number of particles in the initial state, ^ P (e), 
is t he fu nction of energy determined by the first of 
eqs.(|2~T9|). 

For this field configuration to be an (approximate) so- 
lution to the equations of motion, the exponent must be 
stationary with respect to p and T, 



as p = 

dp 
dS p 



(4.2) 



These equations determine the values of p and T in terms 
of energy. 

The field of the periodic instanton is the sum of contri- 
butions of individual constrained instantons. The field of 
a single constrained instanton located at the point t = tq, 
where r = Imi is the Euclidean time, can be approxi- 
mated as follows (throughout this section we use units 
m = 1), 



)(t - t , a;) 



8 Ki(R) 



R 



Here K\ is the modified Bessel function and R = 
\J{ T — T o) 2 + x 2 + p 2 . Hence, the periodic instanton 
field is 

<I>p ( x ) = ^2 ( I ) o(t - nT,x) 



ri, „ sin(fcr)ch(w k r) 
o 



P 



where r = \ x 2 and |r| < T/2. Substituting this expres- 
sion into the Euclidean version of the action fl3.1| ) and 
evaluating the time integral over the period one obtains, 
after some algebra, 

S p {p,T) = ^-^{2\og{p/2) 

+2 7 +l + f(T)}+0(p\pVT 2 ), 



where 



f(T) = 8 J 



k 2 dk 



1 



and 7 = 0.577.. is the Euler constant. 

It is convenient to express all relevant quantities as 
functions of T rather tha n E . The value of p is deter- 
mined by the first of eqs.(4.2), 



p 2 (T)=4exp{-2 7 -2-/(T)}, 



while the second of eqs.(4.2) gives 
32ttV 



E{T) = 



X 



k Z dk- 



l) 2 



(4.4) 



Comparing this to eqs.( 2.19 ) one finds also the number 
of particles, 



N(T) 



327T 2 p 2 

A 



k 2 dk 



u k T 



(4.5) 



The exponential sup pression (as a function of T) is de- 
termined by eq.( |4.1| ). When the period T is expressed 
through energy from eq.(4.4), both the suppression fac- 
tor (4.1) and number of particles (4.5) become functions 
of energy. The latter will be shown in FigM 



V. NUMERICAL RESULTS 

At energies of order sphaleron energy or higher, the 
analytical calculation of the multiparticle rate is not pos- 
sible, and one has to rely on numerical methods. In the 
case of periodic instanton, i.e. at 6 = 0, such calcula- 
tions were performed previously in different models in 
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refs. j29| , [30[ ] . In th e ca se 9^0, the preliminar y re sults in 
the scalar model (3.1) were published in ref. |L9[ . 

The key step of the calculation is the numerical solu- 
tion of the boundary value problem ( 2.18a| ). There are 
several features of this problem which make it difficult 
for numerical analysis. First, eq.(2.18a) is non-linear. 
Second, at 9 ^ th e field is necessarily complex; the 
exponent in eg. ( 2. 17 ) is not positive-definite and the so- 
lution is essentially a saddle point. Third, the contour 
ABC includes both Minkowskian (AB) and Euclidean 
(BC) parts, so the numerical techniques to be used must 
be suitable for both hyperbolic- and elliptic-type equa- 
tions. Finally, the initial boundary conditions (2.18c) 
should be imposed in the region where the field is lin- 
ear, i.e. the physical volume (always finite on the lattice) 
must be large enough in order to allow for linearization. 
The simplification to be mentioned is that eqs.( p.l8"a ) 
are 0(3)-symmetric, so that the problem is actually two- 
dimensional. 

In combination, these factors impose severe constraints 
on both the parameters of the lattice and numerical tech- 
niques to be used. In the lattice version (for details see 
Appendix A) the boundary value problem (2.18a) be- 



smaller than the maximum lattice momentum ~ im x /L, 
which implies 



comes a set of coupled non-linear algebraic equations for 
the values of the field <f>ij = <f>(ti,rj) at the lattice sites 



with coordinates (ti,rj), where ro,... 
the spatial radial direction, while to, 



%x correspond to 
, t nt are complex 
time coordinates lying on the contour ABC in the com- 
plex time plane. The lattice size is characterized by three 
parameters: Tm and T, the Minkowskian and Euclidean 
sizes of the time contour, respectively, and the spatial 
size L. On the Euclidean part of the contour, the so- 
lution is compact in space. On the Minkowskian part, 
when viewed as evolving backward in time from B to A, 
the solution propagates along the light cone |t| = r. The 
spatial size of the lattice is therefore determined not by 
the characteristic size of the Euclidean field configura- 
tion, but by the requirement 

L ;> T M , 



E L 
N 7rn 7 



< 1. 



Thus, n x directly controls the available region of param- 
eter space. 

The numerical method to solve the set of equations 
which constitu te the lattice version of the boundary 
value problem (2.18a) should be chosen according to the 
specifics of the problem described above. To get rid of 
the non-linearity we employ the multidimensional analog 
of the Newton-Raphson method J2jj which approaches 
the desired solution iteratively. At each iteration, the 
linearized equations in the background of the current ap- 
proximation have to be solved. The next approximation 
is obtained by adding the solution to the background, 
and the procedure is repeated. The advantage of the al- 
gorithm is that it does not require positive-definitcness 
of the matrix of second derivatives. It is, however, sen- 
sitive to zero modes. In the absence of zero modes, the 
algorithm converges quadratically; the accuracy of 10 -9 
is typically reached in 3-5 iterations. The convergence 
slows down in the presence of very soft modes, as typi- 
cally happens near bifurcation points. 

The complication is that, i n con tinuous formulation, 
the boundary value problem (2.18a) does have an exact 
zero mode. This zero mode corresponds to translation 
in the real time direction (both field equations and the 
boundary conditions are invariant under such a transla- 
tion). In the lattice version, due to discretization and 
finite volume effects, this zero mode transforms into a 
very soft one. In order to avoid the convergence prob- 
lems related to this soft mode we impose a constraint 
which breaks the translational symmetry (for details see 
Appendix A). 

At each Newton-Raphson iteration one has to solve the 
set of n t x n x linear equations of the general form 



where Tm in turn is determined-by the linearization time. 
The latter is model-dependentE. In the model (3.1) Tm 



can be taken as small as Tm = 2.5m . The correspond- 
ing value of L in our calculations was L = 3m _1 . 

Since the minimal physical volume is fixed in a given 
model, the grid resolution is determined by the number 
of mesh points in time and space directions, n t and n x , 
respectively. In the E-N space, the two regions of inter- 
est are large E at fixed N and small N at fixed E, both 
corresponding to large values of the average momentum 
of initial particles ki n ~ E/N. The latter has to be much 



5 Linearization is reasonably fast in four-dimensional models 
due to the volume factor. In two dimensions it can be very 
slow IE2I1. 



where u is the vector formed of n t n x unknowns, L is 
the matrix of dimension n t n x x ntn x (first variation of 
the full non-linear equations) and d is a constant vector 
(full equations evaluated at the current background; at 
the desired solution d = 0). The matrix L is neither 
positive-definite nor even symmetric, but has a special 
sparse structure as it originates from the second order 
differential equation. The inversion of this matrix is the 
most time consuming part of the calculation; its effective- 
ness determines how large nt and n x can be taken. In 
our calculation we used the forward elimination and back- 
substitution algorithm described in Appendix B, which 
amounts to ~ n t n x multiplications and requires the stor- 
age space for ~ n t n x double precision numbers. Note 
that this algorithm is asymmetric in n t and n x and is 
suitable for the case nt n x . The results presented be- 
low were obtained at n x = 80 and n t = 400 — 500. We 
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have been able to implement a more efficient algorithm 
from the conjugate-gradient family |27f | only in the case 
= 0. 

The Newton-Raphson method requires a good initial 
approximation for the solution. This favors the follow- 
ing general strategy. We first find the periodic instanton 
solution (which corresponds to 9 = 0) with a given pe- 
riod T close to the period of oscillations in the sphaleron 
negative mode. It can be approximated by the sphaleron 
configuration plus oscillation in the negative mode (both 
the sphaleron and its negative mode have to be known 
for a given lattice; the initial amplitude of the oscillation 
has to be chosen by trial and error). After the periodic 
instanton is found, we change parameters T and 9 by 
small steps, using the solution from the previous run as 
a starting configuration for the next one. At each step 
we calculate the energy E, number of particles N and 
the exponential suppression factor F(e,v). By making 
use of smoothness of E and F as functions of T and 9, 
the changes in the parameters can be organized so as to 
maintain E or F, whichever is desired, (approximately) 
constant. 
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FIG. 5. Lines of F{e,v) = const. The value of F changes 
from F = for the uppermost (solid) line to F = 0.6 with 
the step 0.05, in the normalization F(0, 0) = 1. The line 
directed from the sphaleron (e = u = 1) to the zero energy 
instanton (e = v = 0) is formed by the periodic instanton 
solutions. Each point on the dia gram represents the solution 
to the boundary value problem (2.18a). 



The behavi or o f the multiparticle rate we have found 
in the model (3d) is summarized in Fig.^j which shows 
the lines of constant F in the e— v plane. The value of 
F changes from F — for the uppermost (solid) line to 
F = 0.6 for the lowest one, with the step 0.05. Units 
are such that the zero energy instanton suppression is 
F(0,0) = 1. Each dot in Fig|5 repre sents the solution to 
the boundary value problem (2.18a) with corresponding 
e and v. For completeness, Fig. [3] also shows the periodic 
instanton solutions which form the line directed from the 
sphaleron e = v = 1 to the zero energy instanton e = v = 
0. 

The behavior of the function F is rather remarkable. 



Near the periodic instanton, the dependence of F on v 
is very weak. The number of particles in the initial state 
can be substantially reduced with almost no increase in 
the exponential suppression. At the very periodic instan- 
ton, i.e. at 6 = 0, the slope of the lines of constant F be- 
co mes in finite. The latter is easy to see analytically from 
eq. ( 2.17| ) by making use of the fact that F is stationary 
with respect to variations of T and 9. One finds 



dN 
~dE 



The r.h.s. goes to infinity at 9 — ► 0. 

The behavior of the function F changes as one moves 
away from the periodic instanton. At small v or at high 
energies, the dependence of F on energy becomes weak. 
Unlike the vicinity of the periodic instanton, in this re- 
gion the increase of energy does not lead to a noticeable 
reduction in the exponential suppression. 

As follows from Fig|| the dependence of F on the num- 
ber of particles at fixed energy is monotonic, in agree- 
ment with general arguments of ref. |T^,[l4|]. Thus, the 
two-particle cross section 



"/..'If I - r:, V {-^-F(e,0) 



is certainly exponentially suppressed at E 3.5E sp h- It 
is clear from FigJ| that the energy E* (if any) at which 
the exponential suppression may disappear is substan- 
tially higher. A simple estimate is obtained by taking 
the average slope of the line F = in the energy range 
2 < e < 3 and performing linear extrapolation. One finds 
that the extrapolated line F = crosses v = axis at 
the point e* ~ 10. This gives the following lower bound 
for the energy E* , 

E* > 10E sp h- 
This value is likely to be underestimated. 
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FIG. 6. F(e, v) as a function of v at fixed e. Numbers near 
the curves show the values of e. 
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Our present data allow a rough estimate of the expo- 
nential suppression of the two-particle cross section in 
the region E ~ E sp h. This estimate can be obtained by 
extrapolating the function F to v = 0. Fig.|] shows the 
dependence of the function F on is for various values of e 
in the range e ~ 1. Different curves of Fig.^ correspond 
to e = 0.6,0.7,0.8,0.9,1.0,1.1,2.0 from top to bottom. 
The estimated value is 

F HG {e = 1) « 0.8 

with the accuracy of order 10%. Fig.|| confirms the con- 
clusion that the dependence of Fhg{^) on e is slow: all 
curves, including the lowest one corresponding to e = 2, 
converge as v goes to zero and point roughly to the same 
value around 0.8. One concludes that at E ~ E sp h 
the zero energy instanton suppression is reduced by only 
about 20%. For more accurate estimate of the value of 
Fhg(c) more data are needed in the region of small v. 



VI. DISCUSSION AND CONCLUSIONS 

Since no analogous calculations have been performed 
previously, it is worth to discuss in more detail the in- 
terpretation of the numerical results. The first point to 
be considered is their correspondence to the continuum 
theory. The algorithm we used in numerical calculations 
is asymmetric in m and n x ; it mainly restricts the space 
grid resolution. So, the value of n x is our main concern. 




FIG. 7. Low energy perturbation theory for the periodic 
instanton (x) versus numerical results at n x — 80 (o) and 
n x = 240 (solid line) space points. 

The results presented in the previous Section were ob- 
tained on the lattice with n x — 80 space points in the 
space volume L = 3m _1 . In Fig. ^ we plot the line of pe- 
riodic instantons from Fig.|5| (shown by diamonds) versus 
analytical expressions of Sect. 4 (shown by crosses). The 
region of validity of perturbation theory does not over- 
lap with n x = 80 data. Fortunately, in the case 6 = 
the grid resolution can be increased. In Fig]?] the solid 
line represents periodic instantons obtained numerically 



on the lattice with n x = 240 points. This line matches 
both the n x = 80 data and perturbation theory. It is 
clear from Fig|j] that n x = 80 data are accurate above 
v = 0.35, at least in the case 6 = 0. Fig.0 also indi- 
cate s th e region of validity of the perturbative results, 
eqs.@)-§]|). 
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FIG. 8. Lines of constant F at two grid resolutions: 
rix = 80 (o) and n x = 40 (solid lines). 

In the general case 6 ^ we cannot match the per- 
turbative and numerical results because we do not have 
numerical data at sufficiently small v and e. In order 
to check the dependence on grid resolution, we plot in 
Fig.|| the n x = 80 data versus n x = 40 data of ref. [jl9| . 
The n x = 40 data are shown by solid lines. The agree- 
ment is good for v > 0.35 and gets worse for smaller 
v, in accord with our previous estimate. The important 
observation is that the high-resolution curves lie above 
the low-resolution ones. Thus, both the exponential sup- 
pression of the two-particle cross section and the energy 
E* at which it may become exponentially unsuppressed 
are somewhat higher in the continuum limit than follows 
from our lattice results. 

We turn now to the interpretation of the numerical 
solutions presented in Fig.|| as those describing the false 
vacuum decay. As mentioned in Sect. 2, there may be sev- 
eral solutions to the boundary value problem (2.18a); the 
physically relevant one is that continuously connected to 
the low energy periodic instanton. Since Fig.[| demon- 
strates perfectly smooth behavior and suggests the ab- 
sence of bifurcations in the scanned region of parameters, 
we believe that the solutions of Fig.|| indeed describe the 
false vacuum decay. Here we consider more direct argu- 
ment. 

As explained in Sect. 2, the real time part of each solu- 
tion describes the evolution of the system after tunneling. 
By looking at the final field it is possible, in principle, to 
determine whether the false vacuum decay took place. In 
the model ( |3.lD the clear signature of the decay is a sin- 
gularity of the field on the positive part of the real time 
axis (see Figj|). 
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FIG. 9. The position of singularities (shown by crossed cir- 
cles) of the solution to the boundary value problem (2.18a) in 
the —\cj) 4 model. The boundary value problem is solved on 
the contour ABC. 

Since the boundary value problem is solved on the con- 
tour ABC, the real time part of the solution has to be 
found separately by solving (numerically) the initial data 
problem along the real time axis, with initial values of the 
field and its time derivative determined by the solution 
on the contour ABC. We have performed such calcula- 
tions and found that most of the solutions of Fig.|| indeed 
have singularity on the positive part of the real time axis, 
whereas some solutions, in particular those with e <; 1.5, 
apparently do not have such a singularity. There is, how- 
ever, a problem with the interpretation of the results of 
this check. The reason is that the initial data problem 
mentioned above is unstable. 

Let us consider the situation in more detail. As ex- 
plained in Sect. 2, the solutions to the boundary value 



problem (2.18a) have specific analytic structure which 
requires the presence of singularity between the negative 
part of the real time axis and the contour ABC. Iiij-jthe 
model (3.1) this singularity lies at real negative timea as 
shown in Fig.^j. Since the boundary value problem is in- 
variant under translations in the real time direction, only 
the distance between the two singularities of Fig.[| has in- 
variant meaning. In numerical calculations the transla- 
tional invariance is fixed by the constraint in such a way 
as to keep the left singularity at approximately the same 
place far from the asymptotic region A and not too close 
to the origin. All solutions of Fig.^| have left singularity 
in the range —0.5 < t < —0.3. 

As moving from left singularity of Fig.^J to the right 
one along the real time axis, the field comes from infin- 
ity, bounces off the potential barrier and goes back. At 
E < E sp h this always takes finite amount of time. At 
E > E sp h there exist solutions which spend long time os- 



We believe that this is specific to models with the potential 
unbounded from below. 



cillating above the sphaleron. Clearly, these solutions are 
unstable; small perturbations may cause them to roll fi- 
nally to the false vacuum instead of going back to infinity. 
Numerical experiments show that the instability is very 
strong and gets worse for higher energies. Because of this 
instability, the position of the right singularity cannot be 
determined reliably when the distance to the singularity 
is larger than a certain amount, typically ~ 0.7 for low 
energies and ~ 0.5 for high energies. This can be checked 
by studying the dependence of the distance between the 
singularities on the constraint for given E and N . When 
the constraint is imposed in such a way that the singular- 
ities are at the same distance from the origin, their posi- 
tions are determined better. In the asymmetric case, the 
measured distance between singularities is always larger 
and the right singularity may even disappear. 

At low energies and far from the line F = 0, the dis- 
tance between singularities is small enough to be deter- 
mined reliably. In this region the right singularity is al- 
ways present and the distance between the singularities 
does not depend on the constraint. At high energies or 
close to the line F = 0, the distance is larger and be- 
comes constraint-dependent, so that its real value can- 
not be determined by the above method. On the (E — N) 
plane there is no well-defined boundary between the two 
regions. This makes us to infer that the interpretation 
of our solutions as describing the false vacuum decay is 
correct. 

Let us come to conclusions. Numerical calculations we 
have performed demonstrate that the formalism of refs. 
]l3|,pT 15 1 indeed provides a practical way to calculate 
the exponential suppression of the two-particle cross sec- 
tion of induced tunneling, Fhg(^)- The auxiliary object 
which we actually calculated, the logarithm of the multi- 
particle probability F(e,v), has regular behavior consis- 
tent with theoretical expectations. This allows to set 
lower bounds on energy E* at which the two-particle 
cross section may become exponentially unsuppressed, 
and to estimate the value of the suppression at energies 
E ~ E sp h. We believe that analogous calculations are 
possible in any model with tunneling transitions. 

In the — A</> 4 model we have found that the exponential 
suppression of the total cross section of induced false vac- 
uum decay persists at least to energies of order lQE sp h- 
At E ~ E sp h about 80% of the zero energy suppression 
is still present. The accuracy of this estimate crucially 
depends on the grid resolution, as can be seen from FigJ|. 

Finally, let us comment on the similarity between the 
— A^> 4 model and the SU(2) Yang-Mills-Higgs system. It 
arises from the softly broken conformal symmetry present 
in both cases. Due to this (broken) symmetry, the low en- 
ergy tunneling in both cases is dominated by constrained 
instantons whose size tends to zero in the limit E — > 0. 
The multiparticle rates should also have similar behav- 
ior, at least in the low energy domain. Therefore, we 
expect that in the SU(2) Yang-Mills-Higgs case it should 
be possible to calculate the function Fhg(c) numerically 
along the lines we followed in this paper. 
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APPENDIX A: LATTICE FORMULATION OF 
THE BOUNDARY VALUE PROBLEM 

In this Appendix we gi ve the lattice formulation of the 
boundary value problem ( [2.18a ) and describe the method 
which we use to solve it. 

We first note that the boundary value problem ( [2.18a ) 
is 0(3)-symmetric, so we may restrict ourselves to field 



configurations depending on 



x 2 and t. In order to 



avoid the discretization of r -dependent derivative terms 
in the action, it is convenient to change variables accord- 
ing to 

1 



4> = -ip. 

r 



The action becomes 



S = 4?r 



OO 

J dtj dr{^ 2 



(Al) 



In order to obtain the lattice formulation of the boundary 
value problem ( [2.18a ), we first derive the lattice version 
of the exponent in eq.([2.4|) and then repeat all steps of 



Sect. 2 leading to eqs. (2.18a) 



Consider first the discretization of the action ( |A1[ ) . In 
the lattice formulation, it depends on (n t + l)(n x + 1) 
complex variables ipij — ip(ti,rj), where r$ G [0, L] with 
ro = 0, r„ x = L, while U are complex numbers lying on 
the contour ABC of Fig.® so that t = -T M + iT/2, 
t nt = 0. When T is small (this is the case in the high 
energy domain) the contour ABC passes close to the sin- 
gularity of the field and the solution looses the accuracy. 
In this region of parameter space we deform the contour 
(leaving the point t = and the initial asymptotic region 
intact) so as to avoid the singularity. 

Since fields are associated with the lattice sites while 
their derivatives are associated with links, we define two 
different sets of intervals 



dtn 



r j+ i 



3 = 0. 



1, 



and 



drj = {drj^i + dr 3 )/2, j = 1 . . .n x - 1, 



and similarly for dti and dti. The tilted intervals are used 
to approximate integrals with the integrand defined at 
lattice sites, while non-tilted ones are used for integrands 
defined on links. With these definitions, the discretized 
action (|Alf) reads 



3 f> 



(A2) 



where 



This discretization scheme has the advantage of being 
symmetric and additive: the action for the whole lattice 
is the sum of the actions of elementary plaquettes. It has 
the accuracy 0(dr 2 ). The boundaries are treated to the 
same accuracy, which is i mport ant since the derivation 
of the boundary problem (2.18a) requires integrations by 
parts. 



Consider now the boundary term in eq.(2.4). On the 
lattice, the plane waves are no longer exact eigenfunctions 
of the Hamiltonian. To find their analog consider the 
quadratic part of eq. (|A2|) in the limit of continuous time, 



S i2) = 4tt / dt 



b)dr 5 



(V'j+i - V'j) 2 

2du 



\ifidrj 



where ipj = V'(^ r j)- The kinetic term of this action 
takes the canonical form in terms of variables Xj(t) = 

minus free lattice Hamiltonian 
1 x ^ 



The rest of the integrand can be written as 



where 



2 ^ 

jk 



1 



hjkXjXk, 



1 



Ijk 



drjdrj^i drjdrj 
$j+l,k &j-i,k 



1 



dTi \l dTj dr k dr k J drj dr k 



The diagonalization of hj k determines the eigenfunctions 

and eigenvalues tu 2 which substitute the plane waves 
and frequencies uj^. We perform the diagonalization nu- 
merically. In fact, all expressions of Sect. 2 involving mo- 
mentum representation can be readily translated to the 
lattice language by means of the substitutions 
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w k — > w n 



In particular, the boundary term in eq.(2.17) becomes 

Ef 1 - 7 1+7 
1 ~ i — I — Re TO Re Wfc + q Im TO Im ^ofe 
[1+7 1 — 7 

(A3) 



where 




(n) 



The matrix f2j£ is calculated numerically. Note that it 
has to be calculated only once for a given spatial lattice 

Now we are in position to derive the lattice version of 
the boundary value problem ( 2.18a ). The lattice analog 
of the field equation 6S/8<f> = is 



dS 

dipij 







where i = l...nj — 1. The final boundary conditions im- 
mediately translate to 



dS 



Im 



= 0, 

= 0. 



(A5) 



In order to derive the initial boundary conditions, one has 



to consider the lattice version of the exponent in eq.(2.4), 
take the derivative with respect to ipoj and then set ip' j = 
(ipaj)* and 7k = 7 = c~ e as discussed in Sect. 2. The 
result reads 



as ^ f.1-7 



dif)> 



1 + 7 

Re V>ofc - Im V'ofe \ = 



1-7 



(A6) 



Note that total number of equations matches the num- 



ber of unknowns. Eqs.( A4)-( A6) form a set of coupled 
non-linear algebraic equations which consti tute t he lat- 
tice analog of the boundary value problem (2.18a). 

As has been noted in Sect. 5, in continuum version these 
equations are invariant under translation in the real time 
direction, which leads to the existence of almost zero 
mode on the lattice, i.e. co n tinu ous family of (approx- 
imate) solutions to eqs.( A4)-( A6). This spoils the con- 
vergence in the Newton-Raphson method and has to be 
cured. The standard trick is to introduce the constraint 
which breaks the translational invariance. In the vicin- 
ity of the sphalcron it is natural to require that at the 



point t = the field velocity has zero projection on the 
sphalcron negative mode, 




0. 



(A7) 



where 'E^~ J IS the sphaleron negative mode on the lattice. 
Far from the sphalc ron there is no natural choice, so we 
use the constraint (A7) at all values of parameters. In 
the continuum theory physical quantities do not depend 
on the constraint. 



APPENDIX B: NUMERICAL ALGORITHM 

Here we give basic formulae for the forward elimination 
- backward substitution algorithm which we use to solve 
the linear equations arising at each Newton-Raphson it- 
eration. 

The general form of the linearized equations is 



L ■ u = d, 



(Bl) 



where L is the matrix of first derivatives of the full equa- 
tions, Uij are unknowns, while —dij are full equations 
(A4) calculated at the current background (dij = if the back- 
ground is a solution). The matrix L has the dimension 
(n t + l)(n x + 1) x (n t + l)(n x + 1). It is a sparse ma- 
trix with a special structure which is most conveniently 
represented in the following form. Let us introduce vec- 
tor notations in which the spatial index j is implicit, 
Ui = {itio, ■ ■ • ,WmJ and similarly for d{. Then the ma- 
trix L has block-three-diagonal form with the diagonal 



blocks Di and off-diagonal blocks d\ and d[ + \ all 



of the same dimension (n x + 1) x (n x + 1). Eq.(Bl) is 
equivalent to the following set of equations, 



DqUq 

dHu»_i + DiUi 

D< n~ )u n t -1 + D nt U nt 



do, 



D\ +) u t+1 = di, (B2) 



The first and last of these equations represent the bound- 
ary conditions. 

Let us define a set of matrices Ai of dimension (n x + 
1) x (n x + 1) and vectors i — 0. . .n t — 1, by the 
equations 



A t u l+ i + bi 



(B3) 



The first of eqs.(B2) implies 



A = -[£> ]- 1 M +) , 
b = [D^dv. 



Eqs.(B2) and (B3) together give the recursion relation 
for Ai and 6^, 
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>(-) 



b L = [D [ f , A i - 1 +D i ]-' l -[d i -D < r\_ 1 ]. 

At the stage of forward elimination, all Ai and bi are 
calculated and stored (this procedure is equivalent to the 
elimination of the lower block-subdiagonal of the matrix 
L). 



At the last point i = n t , the third of eqs.(B2) and 
eq.(B3) taken at i = rif — 1 determine the vector of un- 
knowns u„ t , 



[D { n ^A nt ^ + D nt ]-'[d nt 



Other unknowns are found by sequential use of eq.( ]B3| ) 
(back-substitution) . 

Clearly, the most time consuming stage is the calcula- 
tion of the matrices Ai. It amounts to matrix inversion 
and subsequent multiplication of the diagonal matrix by 
the result, all repeated n t times. This can be done in 
~ ntn\ complex multiplications p7| with a coefficient 
close to 1. The algorithm requires storage space for n t n^, 
complex numbers. 
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